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Abstract 

We classify all bifurcations from traveling waves to non-trivial time-periodic solutions 
of the Benjamin-Ono equation that are predicted by linearization. We use a spectrally 
accurate numerical continuation method to study several paths of non-trivial solutions 
beyond the realm of linear theory. These paths are found to either re-connect with a 
different traveling wave or to blow up. In the latter case, as the bifurcation parameter 
approaches a critical value, the amplitude of the initial condition grows without bound 
and the period approaches zero. We then prove a theorem that gives the mapping 
from one bifurcation to its counterpart on the other side of the path and exhibits exact 
formulas for the time-periodic solutions on this path. The Fourier coefficients of these 
solutions are power sums of a finite number of particle positions whose elementary 
symmetric functions execute simple orbits (circles or epicycles) in the unit disk of the 
complex plane. We also find examples of interior bifurcations from these paths of already 
non-trivial solutions, but we do not attempt to analyze their analytic structure. 
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1 Introduction 

The Benjamin-Ono equation is a non-local, non-linear dispersive equation intended to de- 
scribe the propagation of internal waves in a deep, stratified fluid [H [151 EO]- In spite of 
non-locality, it is an integrable Hamiltonian system with meromorphic particle solutions 
[12\ I13j. A^-soliton solutions [24J, and A^-phase multi-periodic solutions [32} 116 1 [26], A bilin- 
ear formalism [32] and a Backlund transformation [28[ [TJ [25] have been found to generate 
special solutions of the equation, and, in the non-periodic setting of rapidly decaying initial 
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conditions, an inverse scattering transform has been developed \18\ [20] that exploits an in- 
teresting Lax pair structure in which the solution plays the role of a compatibility condition 
in a Riemann-Hilbert problem. 

It is common practice in numerical analysis to test a numerical method using a problem 
for which exact solutions can be found. Our initial interest in Benjamin-Ono was to serve 
as such a test problem. Although many of the tools mentioned above can be used to study 
time-periodic solutions, they do not generalize to problems such as the vortex sheet with 
surface tension [21 H] or the true water wave [311 |19j, which are not known to be integrable. 
Our goal in this paper is to develop tools that will generalize to these harder problems 
and use them to study bifurcation and global reconnection in the space of time-periodic 
solutions of B-0. Specifically, we employ a variant of the numerical continuation method 
we introduced in [3] for this purpose, which yields solutions that are accurate enough that 
we are able to recognize their analytic form. 

Because we approached the problem from a completely different viewpoint, our de- 
scription of these exact solutions is very different from previously known representations 
of multi-periodic solutions. Rather than solve a system of non-linear algebraic equations 
at each x to find u{x,t) as was done in [26] , we represent u{x,t) in terms of its Fourier 
coefficients Cfc(t), which turn out to be power sums Ck = 2[(3i + ■ ■ ■ + of a collection of 
N particles f3j{t) evolving in the unit disk of the complex plane as the zeros of a polyno- 
mial z I— >• P{z,t) whose coefficients execute simple orbits (circles or epicycles in C). The 
connection between the new representation and previous representations will be explored 
elsewhere [37J. 

Many of our findings on the structure of bifurcations and reconnections in the manifold 
of time-periodic solutions of the Benjamin-Ono equation are likely to hold for other systems 
as well. One interesting pitfall we have identified by applying our method to an integrable 
problem is that degenerate bifurcations can exist that are not predicted by counting linearly 
independent, periodic solutions of the linearization about traveling waves. Although it is 
possible that such degeneracy is a consequence of the symmetries that make this problem 
integrable, it is also possible that other problems such as the water wave will also possess 
degenerate bifurcations that are invisible to a linearized analysis. We have also found that 
one cannot achieve a complete understanding of these manifolds of time-periodic solutions by 
holding e.g. the mean constant and varying only one parameter. In some of the simulations 
where we hold the mean fixed, the solution (i.e. the norm of the initial condition) blows up 
as the parameter approaches a critical value rather than reconnecting with another traveling 
wave. However, if the mean is simultaneously varied, it is always possible to reconnect. 
Thus, although numerical continuation with more than one parameter is difficult, it will 
likely be necessary to explore multi-dimensional parameter spaces to achieve a thorough 
understanding of time-periodic solutions of other problems. 

On the numerical side, we believe our use of certain Fourier modes of the initial condi- 
tions as bifurcation parameters will prove useful in many other problems beyond Benjamin- 
Ono. We also wish to advocate the use of variational calculus and optimal control for 
the purpose of finding time-periodic solutions (or solving other two point boundary value 
problems). For ODE, a competing method known as orthogonal collocation (e.g. as imple- 
mented in AUTO [17]) has proved to be a very powerful technique for solving boundary 
value problems. This approach becomes quite expensive when the dimension of the system 
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increases, and is therefore less competitive for PDE than it is for ODE. For PDE, many 
authors do not attempt to find exact periodic solutions, and instead point out that typical 
solutions of certain equations do tend to pass near their initial states at a later time [11] , 
If true periodic solutions are sought, a more common approach has been to either iterate 
on a Poincare map and use stability of the orbit to find time-periodic solutions [TD], or use 
a shooting method [331 El] to find a fixed point of the Poincare map. 

In a shooting method, we define a functional F{uq,T) = [u{-,T) — uq] that maps initial 
conditions and a supposed period to the deviation from periodicity. The equation F = 
is then solved by Newton's method, where the Jacobian J = DF is either computed using 
finite differences |35j or by solving the variational equation repeatedly to compute each 
column of J. We have found that it is much more efficient (by a factor of the number of 
columns of J) to instead minimize the scalar functional G = ^H-Fp via a quasi-Newton 
method in which the gradient DG is computed by solving an adjoint PDE. 

Bristeau et. al. [8] developed a similar approach for linear (but two- or three-dimensional) 
scattering problems. Three dimensional problems are intractable by the standard shooting 
approach as J could easily have 10^ columns. However, the gradient of G can be computed 
by solving a single adjoint PDE. The success of the method then boils down to a question 
of the number of iterations required for the minimization algorithm to converge. For linear 
problems, Bristeau et. al. have had success using conjugate gradients to minimize G. We 
find that BEGS [9j works very well for nonlinear problems like the Benjamin-Ono equation 
and the vortex sheet with surface tension [3]. 

To find non-trivial time-periodic solutions in the present work, we use a symmetric 
variant of the algorithm described in [3]. Although the original method works well, we use 
the symmetric variant for the simulations in this paper because evolving to T/2 requires 
half the time-steps and yields more accurate answers (as there is less time for numerical 
roundoff error to corrupt the calculation). Moreover, the number of degrees of freedom in 
the search space of initial conditions is also cut in half and the condition number of the 
problem improves when we eliminate phase shift degrees of freedom via symmetry rather 
than including them in the penalty function described in Section 13.11 Although we do 
not make use of it, there is a procedure known as the Meyer-Marsden-Weinstein reduction 
[271 [23] that allows one to reduce the dimension of a symplectic manifold on which a group 
acts symplectically. This allows one to eliminate actions of the group (e.g. translations) 
from the phase space. Equilibria and periodic solutions of the reduced Hamiltonian system 
correspond to (families of) relative equilibria and relative periodic solutions [39] of the 
original system. 

This paper is organized as follows. In Section[2l we discuss stationary, traveling and par- 
ticle solutions of B-0, linearize about traveling waves, and classify all bifurcations predicted 
by linear theory from traveling waves to non-trivial time-periodic solutions. Some of the 
more technical material from this section is given in Appendix |Al In Section O we present 
a collection of numerical experiments using our continuation method to follow several paths 
of non-trivial solutions beyond the realm of linear theory in order to formulate a theorem 
that gives the global mapping from one traveling wave bifurcation to its counterpart on 
the other side of the path. In Section HI we study the behavior of the Fourier modes of 
the time-periodic solutions found in Section [3] and state a theorem about the exact form 
of these solutions, which is proved in Appendix [Bj Finally, in Section \5\ we discuss inte- 
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rior bifurcations from these paths of already non-trivial solutions to still more complicated 
solutions. Although the existence of such a hierarchy of solutions was already known [32], 
bifurcation between various levels of the hierarchy has not previously been discussed. 

2 Bifurcation from Traveling Waves 

In this section, we study the linearization of the Benjamin-Ono equation about stationary 
solutions and traveling waves by solving an infinite dimensional eigenvalue problem in closed 
form. Each eigenvector corresponds to a time-periodic solution of the linearized equation. 
The traveling case is reduced to the stationary case by requiring that the period of the 
perturbation (with a suitable spatial phase shift) coincide with the period of the traveling 
wave. The main goal of this section is to devise a classification scheme of the bifurcations 
from traveling waves so that in later sections we can describe which (local) bifurcations are 
connected together by a global path of non-trivial time-periodic solutions. 

2.1 Stationary, Traveling and Particle Solutions 

We consider the Benjamin-Ono equation on the periodic interval M/27rZ, namely 

Ut = HUxx - UUx. (1) 

Here H is the Hilbert transform, which has the symbol H{k) = — isgn(A;). The Benjamin- 
Ono equation possesses solutions [121 E] of the form 

TV 

u(x,t) = ao + ^<^(x;A(t)), (2) 
1=1 

where oq is the mean, /3i(t), . . . , ^Nit) are the trajectories of N particles evolving in the 
unit disk A of the complex plane and governed by the ODE 

/^'= E ^^+E^^+^(2^-l-"o)A, {1<1<NI (3) 

m=l P"" m=l PI ~ ^™ 

and (j){x;l3) is the function with Fourier representation 

C 0, A; = 1 

4){k-p) = \ 2/3^ A;>0 i, /3 G A = {z : |z| < 1}. (4) 
y 2p\^\ A; < J 

The function (j){x] f3) has a peak centered at x = arg(/3) with amplitude growing to infinity 
as 1/3 1 approaches 1. The A^-hump traveling waves (with a spatial period of 2tt/N) are a 
special case of the particle solutions given by ([2]) and 

N 

Utr.Ax,t;ao,N,P)=ao + ^cj){x;(3i{t)), /3,(t) = V^e"^^*, c = ao - iVa(/3). (5) 
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Each j3i is assigned a distinct A^th root of /3 and a(/3) is the mean of the one-hump stationary 
solution, namely 

1-3|/3|2 l-a(/3) 



a(/3) 



1-|/3|2' 3-a(/3)- 
The solution (j5]) moves to the right when c > 0. Indeed, it may also be written 

Utray{x,t;aQ,N,l3) = Usta.t{x - ct;N,f3) + c, 

where Ustat is the A^-hump stationary solution 

u,t,t{x;N,/3)=Na{(3)+ H^'^^) = Na{f3) + N<p{Nx; 



(6) 



(7) 



(8) 



The Fourier representation of ttgtat is 



Na{P), k = 0, 

2Nj3''/^, k £ NZ, k>0, 

2NB\^\/^, k £ NZ, k<0, 

otherwise. 



(9) 



Amick and Toland have shown [5] that all traveling waves of the Benjamin-Ono equation 
have the form ([7]); see also [37] . 

2.2 Linearization about Stationary Solutions 

Let u{x) = Ustatix; N, 13) be an A^-hump stationary solution. In we solved the lineariza- 
tion of ([!]) about n, namely 

vt = Hvxx — {uv)x = iBAv, A = Hdx — u, B = —dx, (10) 

i 

by substituting v{x,t) = Re{C2;(x)e*'^*} into ([10]) and solving the eigenvalue problem 

BAz = uJZ (11) 
in closed form. Specifically, we showed that the eigenvalues u}N,n are given by 

500 r 



<^N,n = < 



—^N-n, 




n < ' 









n = 




' (n)(A^- 


n), 


1 < n < A^ - 1 




. (n + 1- 


N){n + 1 + N{1- 


am), n>N ^ 


0, 



N=20, p=1/2 . 



(12) 
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The zero eigenvalue u^Nfi = has geometric multiplicity two and algebraic multiplicity 
three. The eigenfunctions in the kernel of BA are 



,(1,0) 



0^(2;) = -^^istat(a;;A^,/3), ^/^,o(^) 



d 



(13) 
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which correspond to changmg the phase or amphtude of /3 in the underlying stationary 
solution. There is also a Jordan chain [50] of length two associated with z^^'q\x), namely 



,(1,1), 



x) 



iBAz 



(1,1) 

N,0 



J1,0) 
^N,0 



(14) 



which corresponds to the fact that adding a constant to a stationary solution causes it to 
travel. The fact that all the eigenvalues ico^^n in the linearization (jlOp are purely imaginary 
is a consequence of the Hamiltonian structure [13] of the Benjamin-Ono equation. For non- 
Hamiltonian systems, one does not generally expect to find time-periodic perturbations of 
traveling waves (as periodic solutions of the linearized problem may not even exist). 

The eigenfunctions z^^n{x) corresponding to positive eigenvalues ujj\i^n (with n > 1) have 
the Fourier representation 



%n(fc) 



k=n+jN 



C 1 + 



Nj 



/SJ+i j>0 



c 



1 < n < iV- 1 \ 

. -nN 

' (Af-n) [n+{N-n)\p\^\ j ' 



ZN,n{k) 



k=n+l~N+jN 



1 



1 



N 

71+1 



1 



n+l 




(n > N), (15) 



with all other Fourier coefficients equal to zero. The eigenfunctions corresponding to nega- 
tive eigenvalues uJN,n (with n < —1) satisfy 2;Ar^„(x) 
appear in reverse order, conjugated. For 1 < n < 



ZAr^_n(x), so the Fourier coefficients 
1, any linear combination of z^^nix) 
and ZN,N-n{x) is also an eigenfunction; however, the choices here seem most natural as 
they simultaneously diagonalize the shift operator (discussed below) and yield directions 
along which non-trivial solutions exist beyond the linearization. Said differently, we have 
listed the first — 1 positive eigenvalues ujN,n in an unusual order (rather than enumerat- 
ing them monotonically and coalescing multiple eigenvalues) because this is the order that 
leads to the simplest description of the global paths of non-trivial solutions connecting these 
traveling waves. 



2.3 Classification of bifurcations from traveling waves 



Time-periodic solutions of the Benjamin-Ono equation with period T have initial conditions 
that satisfy F{uq,T) = 0, where F : xR ^ is given by 



F{uo,T) = u{-,T) - uq, ut = Huxx-uUx, u{-,0)=uo. 



(16) 



We begin by linearizing F about an A^-hump stationary soution uo{x) = Usta.t{x] N, /?). The 
Frechet derivative DF = {DiF, D2F) : x M — )• yields directional derivatives 



DiF{uQ,T)vo 



d_ 

de 
d 



£=0 



F(no + e^o, T) = v{;T) - vo = [e*^^^ - /] vo, 



(17) 



D2F{uo,T)t = J- F{uo,T + et) 



0. 



6 



Note that vq £ ]ier DiF{u,T) if and only if the solution v{x,t) of the linearized problem is 
periodic with period T. As a result, a basis for the kernel M = ker DF{uo,T) consists of 
(0; 1) together with all pairs {vq;0) of the form 

vo{x) = Re{zN,n{x)} or voix) = lm{zN,n{x)}, (18) 

where n ranges over all integers such that 

u;N,nT e 2ttZ (19) 

with and /3 (in the formula (jl2p for iON,n) held fixed. The corresponding periodic solutions 
of the linearized problem are 

v{x,t) =Re{zN,nix)e''^'''"^} or t>(x, t) = Im{zAr,„(2;)e*'^^-"*}. (20) 

Negative values of n have already been accounted for in (fT8|) and ([20|) using zj^^^n{x) = 
ZN,n{x), and the n = case always yields two vectors in the kernel, namely those in (I13p . 
These directions do not cause bifurcations as they lead to other stationary solutions. 

Next we wish to linearize F about an arbitrary traveling wave. Suppose u{x) = 
Usta.t{x; N, 13) is an iV-hump stationary solution and U{x,t) = u{x — ct) + c is a travel- 
ing wave. Then the solutions v and V of the linearizations about u and U, respectively, 
satisfy V{x,t) = v{x — ct,t). Note also that 

F{Uo,T) = iff cT='^ for some i/eZ, (21) 

where Uo{x) = U{x,0) = u{x) + c. Note that v is the number of times the traveling wave 
turns over itself in one period. Assuming (j2ip holds, we set 6 = 2'ku/N and compute 



(22) 



[DiF{U^,T)vo\{x) = v{x - cT,T) - vo{x) = [(^^e^^^^ - I)vo]{x), 
[D2F{Uo,T)t]{x) = Utix,T)T = -cux{x - cT)t = -cUx{x)t, 

where v solves (jlOp and the shift operator Sg is defined via 

Sez{x) = z{x - e), Se,ki = e-^^'Ski- (23) 

One element of = ker DF(Uo,T) arises from ()14p . which gives 

e'^^*l = l-iux ^ DiFiUo,T){-c/T) + D2F{Uo,T)l = 0, 

and implies {—c/T; 1) G N . This just means that we can change the period T by a small 
amount r by adding the constant — (c/T)r to Uq; (this also follows from the condition (j2ip 
that cT = 6 = const). If we wish to change the period without changing the mean, we need 
to simultaneously adjust |/3| in the underlying stationary solution u{x) = Usta.t{x; N, f3). The 
other elements of M are of the form {vq; 0) with 

vo{x) = Re{zN,n{x)} or vo{x) = lm{zN,n{x)}. (24) 

The admissible values of n here are found using (I22p together with 

S,e^^^^.^,. = e^(---^-^'=--)z;v,n, e = (25) 
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where k]\f^n is the stride offset of the non-zero Fourier coefficients of z^^n, i-e. 

ZN,nik) / k- kN,n £ NZ. (26) 

Thus, instead of (|19p . n ranges over all integers such that 



-kN-n, n < 0, 

n = 0, 

n 1 < n < iV 

mod(n+l,A^) n>N. 



(27) 



As before, negative values of n need not be considered once we take real and imaginary 
parts in (plj) . and the n = case always gives the two vectors (z^'q^; 0), (z^q; 0) € AA, 
which lead to other traveling waves rather than bifurcations to non-trivial solutions. 

Our numerical experiments have led us to the following conjecture, which we prove as 
part of Theorem [3] in Section |4l 

Conjecture 1 For every /3 € A and {N, u, n, m) € satisfying 

N >l, u eZ, n > 1, m > 1, m G ukN,n + NZ, (28) 

there is a four parameter sheet of non-trivial time-periodic solutions bifurcating from the N- 
hump traveling wave with speed index u, (cT = 2'ku/N), bifurcation index n, and oscillation 
index m, (uj^^nT = ^ixrajN). The phase and amplitude of the traveling wave are determined 
by p. 

The main content of this conjecture is that we do not have to consider linear combi- 
nations of the ZM,n with different values of n to find periodic solutions of the non-linear 
problem — this basis is already "diagonal" with respect to these bifurcations. This is true 
in spite of a small divisor problem preventing DF(JJo, T) from being Fredholm. The decision 
to number the first A^ — 1 eigenvalues coN^n non-monotonically in (jl2p and to simultaneously 
diagonalize the shift operator Sg when choosing eigenvectors ZN,n in (flSl) was essential to 
make this work. Formulas relating the period, T, the mean, ao, and the decay parameter, 
1/3 1, for each of these bifurcations are given in Appendix [Al along with a list of bifurcation 
rules governing "legal" values of the mean. 

A canonical way to generate one of these bifurcations is to take /3 real and perturb the 
initial condition in the direction vo{x) = 'Re{zN,n{x)}. This leads to non-trivial solutions 
with even symmetry at t = 0. Perturbation in the lm{z]\f^n{x)} direction yields the same 
set of non-trivial solutions, but with a spatial and temporal phase shift: 

lm{zNA^ - ct)e*'^*} = Re {zN,n ((^ " " " i;)) } ' (29) 

where cj = ujn^u- The manifold of non-trivial solutions is four dimensional with two essential 
parameters (e.g. the mean ao and a parameter governing the distance from the traveling 
wave) and two inessential parameters (the spatial and temporal phase). In our numerical 
studies, we use the real part of a Fourier coefficient Ck of the initial condition (with k such 
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that ZN,n{k) 7^ 0) for the second essential bifurcation parameter. When we discuss exact 
solutions in Section HI a different parameter will be used. 

We remark that this enumeration of bifurcations accounts for all time-periodic solutions 
of the linearization about traveling waves; therefore, the heuristic that each bifurcation 
of the non-linear problem gives rise to a linearly independent vector in the kernel M of 
the linearized problem suggests that we have found all bifurcations from traveling waves. 
Interestingly, this turns out not to be the case; the interior bifurcations we discuss in 
Section[5]can occur at the endpoints of the path, allowing for degenerate bifurcations directly 
from traveling waves to higher levels in the infinite hierarchy of time-periodic solutions. Only 
the transition from the first level of the hierarchy to the second is "visible" to a linearized 
analysis about traveling waves. The other transitions become linearly dependent on these 
in the limit as the traveling wave is approached; they will be analyzed in [37] . 

3 Numerical Experiments 

In this section we present a collection of numerical experiments in which we start with a 
given bifurcation (N, u, n, m, (3) and use a symmetric variant of the method we described 
in [3j for finding periodic solutions of non-linear PDE to continue these solutions until 
another traveling wave is found, or until the solution blows up as the bifurcation parameter 
approaches a critical value. We determine the bifurcation indices {N\ v', n', ml) at the other 
end of the path of non-trivial solutions by fitting the data to the formulas of the previous 
section. By trial and error, we are then able to guess a formula relating (N' ,1^' ,n' ,m') to 
{N ,v^n,m) that we use in Section H] to construct exact solutions. 

3.1 Numerical Method 

As mentioned in Section 12.31 a natural choice of spatial and temporal phase can be achieved 
by choosing the parameter (3 of the traveling wave to be real and perturbing the initial 
condition in the direction vq{x) = Re{2;Ar^„(x)}. For reasons of efficiency and accuracy 
(explained in the introduction), we now restrict our search for time-periodic solutions of ([T]) 
to functions u{x, t) that possess even spatial symmetry at t = 0. If we succeed in finding 
solutions with this symmetry, then they (together with their phase-shifted counterparts 
analogous to (j29p ) span the nullspace N = ker DF{Uq,T) in the limit that the perturbation 
goes to zero. Thus, we do not expect symmetry breaking bifurcations from traveling waves 
that cannot be phase shifted to have even symmetry at t = 0. 

The Benjamin-Ono equation has the property that if u{x,t) is a solution of ([T]), then 
so is U{x,t) = u{—x, —t). As a result, if n is a solution such that u{x,T/2) = U{x, — T/2), 
then u{x,T) = U{x,0), i.e. u is time-periodic if the initial condition has even symmetry. 
Thus, we seek initial conditions uq with even symmetry and a period T to minimize the 
functional 



Gtot{uO,T) = G{uo,T) + Gpena.lty{uo,T), 



(30) 



where 




^dx 



(31) 



9 



and Gpcna.ity{uo,T) is a non-negative penalty function to impose the mean and set the 
bifurcation parameter. To compute the gradient of G with respect to variation of the initial 
conditions, we use 



d 
de 



where the variational derivative 



G{uo + evo,T) = [ 

^0 



27r 

- — {x)vo{x)dx, (32) 
duo 



6G 

— (x) = 2w(x,T/2), wo{x) = u{x,T/2) -u{2tt - x,T/2) (33) 
duo 

is found by solving the following adjoint equation from s = to s = T/2: 

Ws{x,s) = -Hwxx{x,s) + u{x,'^ - s)wx{x,s), w{-,^)=wo. (34) 
Since vo is assumed symmetric in this formulation, (|33p is equivalent to 



SG 

— {x) = w{x,T/2) + w{27r-x,T/2). (35) 
duo 

The Benjamin-Ono and adjoint equations are solved using a pseudo-spectral collocation 
method employing a fourth order semi-implicit additive Runge-Kutta method [Ml [211 El] 
to advance the solution in time. The BFGS method [9l [29] is then used to minimize Gtot 
(varying the period and the Fourier coefficients of the initial conditions). We use the penalty 
function 

Gpenalty(uO,T) = ^([ao(0) - Oq]' + [ai^(O) - pf^ (36) 

to specify the mean ao and the real part p of the Kih Fourier coefficient of the initial 
condition 

M/2 

uo{x)= Cfc(0)e^'=^ Ck{t) = ak{t) + ibk{t). (37) 

k=-M/2+l 

The parameters qq and p serve as the bifurcation parameters while the phases are deter- 
mined by requiring that the solution have even symmetry at t = 0. We generally choose K 
to be the first k > 1 such that zj\f^ri{k) 7^ 0. 

Our continuation method consists of three stages. First, we choose a traveling wave 
and a set of bifurcation indices to begin the path of non-trivial solutions. We also choose 
a direction in which to vary the bifurcation parameter p and the mean ao- In most of 
our numerical experiments, we hold ao fixed; however, in the example of Figure [U] below, 
we vary p and ao simultaneously. The traveling wave serves as the zeroth point on the 
path. The initial guess for the first point on the path is obtained by perturbing the initial 
condition of the traveling wave in the direction Re{zAr^„(a;)}. We use the period T given in 
(]77p in Appendix[Slas a starting guess. We then use the minimization algorithm to descend 
from the starting guess predicted by linear theory to an actual time-periodic solution. The 
second stage of the continuation algorithm consists of varying p (and possibly ao)) using 
linear extrapolation for the starting guess (for uq and T) of the next solution, and then 
minimizing Gtot to find an actual time-periodic solution with these values of p and ao- If 
the initial value of Gtot from the extrapolation step is too large, we discard the step and 
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try again with a smaller change in p and oq- The final stage of the algorithm consists of 
identifying the reconnection on the other side of the path. We do this by blindly overshooting 
the target values of p and (which we do not know in advance). Invariably, the algorithm 
will lock onto a family of traveling waves once we reach the end of the path of non-trivial 
solutions. We look at the Fourier coefficients of the last non-trivial solution before the 
traveling waves are reached and match them with the formulas for Z7v',n'(^) to determine 
the correct bifurcation indices on this side of the path. (A prime indicates indices for the 
bifurcation at the other end of the path.) We then recompute the last several solutions on 
the path of non-trivial solutions with appropriate values of p and to arrive exactly at 
the traveling wave on the last iteration. We sometimes change K in ()36p to compute this 
reconnection to avoid Z]\f'^n'{K) = 0. 

The running time of our algorithm (on a 2.4 GHz desktop machine) varies from a few 
hours to compute one of the paths labeled a-l in ([38]) ~ (|1T]) below, to a few days to compute 
a path in which the solution blows up, such as the one shown in Figure [5] below. We 
always refine the mesh and timestep enough so that the solutions are essentially exact 
(with Gtot < 10~^^ in the easy cases and 10"^*^ in the hard cases). 

3.2 Global paths of non-trivial solutions 

We now investigate the global behavior of non-trivial solutions that bifurcate from arbi- 
trary stationary or traveling waves. We find that these non-trivial solutions act as rungs 
in a ladder, connecting stationary and traveling solutions with different speeds and wave- 
lengths by creating or annihilating oscillatory humps that grow or shrink in amplitude until 
they become part of the stationary or traveling wave on the other side of the rung. In 
some cases, rather than re-connecting with another traveling wave, the solution blows up 
(i.e. the L-^-norm of the initial condition grows without bound) as the bifurcation parameter 
p approaches a critical value. However, even in these cases a re-connection with another 
traveling wave does occur if, in addition to p, we vary the mean, oq, in an appropriate way. 

Recall from Section 12.31 that we can enumerate all such bifurcations by specifying a 
complex parameter /3 in the unit disk A along with four integers {N, i/, n, m) satisfying (I28p . 
and in most cases we can solve for |/3| in terms of the mean, ao, using ()80p in Appendix [Al 
In [3], we presented a detailed study of the solutions on the path connecting a one- hump 
stationary solution to a two- hump traveling wave moving left. We denote this path by 

a: (1,0,1,1) ^ (2,-1,1,1), (38) 

where the label a refers to the bifurcation diagram in Figured! We have also computed the 
next several bifurcations {n = 2, 3, 4) from the one-hump stationary solution and found that 
they connect up with a traveling wave with A^' = n -|- 1 humps moving left with speed index 
u' = —1, where we denote the bifurcation on the other side of the path by {N' ,1^' ,n' ,m'). 
By comparing the Fourier coefficients of the last few non-trivial solutions on these paths 
to those of the linearization about the A^'-hump traveling wave, we determined that the 
bifurcation and oscillation indices satisfy n' = n and m' = 1, respectively. Studying these 
reconnections revealed that the correct way to number the eigenvalues lo^i ^n' was to split the 
double eigenvalues with n' < N' apart as we did in (jl2p by simultaneously diagonalizing the 
shift operator and ordering the 0JN',n' via the stride offset of the corresponding eigenvectors 
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First several bifurcations connecting stationary and traveling waves 




T 

Figure 1: Paths of non-trivial solutions listed in equations (|38 p - (j4ip . The second Fourier 
mode of the eigenvector ZN^n{x) in the linearization is non-zero for the pitchfork bifurcations 
and is zero for the one-sided, oblique-angle bifurcations. The point labeled P corresponds 
to the solution in Fig. [3] below. 

(rather than monotonically) . Using this ordering, the non-trivial solutions connect up with 
the A^'-hump traveling wave along the z^'^n' direction (without involving ZN',N'-n')- These 
results are summarized as 

b: (1,0,2,1) ^ (3,-1,2,1), 

c: (1,0,3,1) ^ (4,-1,3,1), (39) 
d: (1,0,4,1) ^ (5,-1,4,1). 

The labels a, b, c, d in ([38|) and (p9|) correspond to the paths labeled 7d, 8d, 5c, a, etc. in the 
bifurcation diagram. When an integer p precedes a label, it means that the period T that 
is plotted is p times larger than the fundamental period of the solution represented. Thus, 
curve 7d is the image of curve d (not shown) under the linear transformation (T, 02) 1— )• 
(7r, 02). In our labeling scheme, we just need to multiply v, m, u', m' by p to obtain the 
new path, e.g. 

7d: (1,0,4,7) ^ (5,-7,4,7). (40) 

In this diagram, we plot 02(0) vs. T with the spatial and temporal phases chosen so the 
solution is even at t = 0. For example, on path d, as we decrease p = 02(0) from 0.371087 
to 0, the solution transitions from the one- hump stationary solution to the five- hump left- 
traveling wave as shown in Figure El 
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Figure 2: Periodic solution on path d connecting the one-hump stationary solution to the 
five-hump left-traveling wave {oq = 0.544375). The second Fourier mode of 2:1,4(3;) is zero, 
which explains why 02(0) = 0.366113 for this solution is only 1.35% of the way between the 
stationary solution 02(0) = 0.371087 and the five-hump traveling wave 02(0) = 0. 



It is interesting that the paths labeled a and 3b in Figure [J meet the one- hump stationary 
solutions in a pitchfork while the other paths (such as 5c and 8d) meet at an oblique angle 
from one side only. This is because the second Fourier mode of the eigenvector zi,n(x) in 
the linearization about the stationary solution is zero in these latter cases, so the change in 
02(0) from that of the stationary solution (namely 0.371087) is a higher order effect, (as is 
the change in T). This explains the oblique angle. We now explain why these bifurcations 
occur from one side only. When we go beyond the linearization as we have here, we find that 
C2{t) = a2{t) +ib2{t) has a nearly circular (epitrochoidal) orbit in case a, a circular orbit in 
case 6, and remains constant in time in cases c and d (see Section |3| . If one branch of the 
pitchfork corresponds to 02(0), the other is a2(T/2) since the function u{-,T/2) also has 
even symmetry. But in cases c and d, 02(0) is equal to a2(T/2) even though the functions 
it(-,0) and u{-,T/2) are different. These cases also become pitchforks when a different 
Fourier coefficient ai^(O) is used as the bifurcation parameter. 

Next we compute the first several bifurcations from the two-hump traveling waves with 
mean ao = 0.544375 and speed index u = —1. We set A'" = 2, = —1, n e {1,2,3,4} 
and choose the first several legal m values, i.e. values of m that satisfy the bifurcation 
rules of Figure [T4l in AppendixO For example, the curves labeled i, j, k and / in Figured 
correspond to the bifurcations (2, —1,4, m) with m = 11, 13, 15, 17; smaller values (and even 
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values) of m are not allowed. In addition to the path a in (j38p above, we obtain the paths 



e : 


(2, 


-1,2,3) ^ 


(3, 


-3,1,3), 


i : 


(2, 


-1,4,11) ^ 


(5, 


-8,3,11), 


/: 


(2, 


-1,3,6) ^ 


(4, 


-5,2,6), 


j ■ 


(2, 


-1,4,13) ^ 


(5, 


-9,3,13), 


g ■ 


(2, 


1,3,8) ^ 


(4, 


-6,2,8), 


k : 


(2, 


-1,4,15) ^ 


(5, 


-10,3,15), 


h : 


(2,- 


1,3,10) ^ 


(4, 


-7,2,10), 


I : 


(2, 


-1,4,17) ^ 


(5, 


-11,3,17). 



(41) 



The paths /, g and h meet the curve representing the two- hump traveling waves in a 
pitchfork bifurcation while the others meet obliquely from one side. This, again, is an 
anomaly of having chosen the second Fourier mode for the bifurcation parameter. The 
dotted line near the path e is the curve obtained when e is reflected across the T-axis. 
Solutions on this dotted line correspond to solutions on path e shifted by 7r/2 in space, 
which changes the sign of p = 02(0) but also breaks the even symmetry of the solution at 
t = 0. The paths labeled i, j, k and / are exactly symmetric when reflected about the T-axis 
because C2{t) has a circular orbit centered at zero in these cases. It is interesting that so 
many of the paths in this bifurcation diagram terminate when T = vr (or a simple rational 
multiple of vr). This is due to the fact that T in ([77|) in Appendix Rl is independent of a 
when n < N. 

The solutions u{x, t) corresponding to points along the paths 6, c and d are qualitatively 
similar to each other. As shown in Figure O these solutions look like A^'-hump waves trav- 
eling over a stationary one-hump carrier signal. At one end of the path the high frequency 
wave may be viewed as a perturbation of the one-hump stationary solution, while at the 
other end of the path it is more appropriate to regard the stationary solution as the pertur- 
bation, causing the traveling wave to bulge upward as it passes near x = tt and downward 
near x = and x = 2tt. In all these cases, the solution repeats itself when one of the high 
frequency waves has moved left one slot to assume the shape of its left neighbor at t = 0. 

By contrast, the solutions that bifurcate from the two-hump traveling waves, i.e. those 
on the paths listed in (j4ip . have the property that when a wave has moved left one slot to 
the location that its neighbor occupied at t = 0, it has acquired a different shape and must 
keep progressing a number of slots before it finally lines up with one of the initial waves. 
This is illustrated in Figure [3] for the solution labeled P in Figure [T] on the path 

e: (2,-1,2,3) ^ (3,-3,1,3). (42) 

This solution is qualitatively similar to the linearized solution (3, —3, 1,3). There are A' = 3 
humps oscillating with the same amplitude but with different phases as they travel left. 
They do not line up with the initial condition again until they have traveled three slots 
(v' = —3) and progressed through one cycle {ml /N' = 3/3), which leads to a braided effect 
when the time history of the solution is plotted on one graph. All the solutions on path e are 
irreducible in the sense that there is no smaller time T in which they are periodic (unlike the 
cases labeled 36, 5c, Id etc. in Figuredl which are reducible to 6, c and d, respectively). Note 
that although u' = —3 and m' = 3 are both divisible by 3, we cannot reduce (3, —3, 1, 3) to 
(3, —1, 1, 1) as the latter indices violate the bifurcation rules of Figure [TH in Appendix lAl 
We also mention that at the beginning of the path, near (2,-1,2,3), the braiding effect 
is not present; instead, the solution can be described as two humps bouncing out of phase 
as they travel left. In one period, they each travel left one slot {v = —1) and bounce 1.5 
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Periodic solution between two and three hump traveling waves 




Figure 3: Time-periodic solution (labeled P in Fig. [T]) on path e connecting two- and three- 
hump traveling waves. The amplitude of each hump oscillates as it travels left. The dotted 
curves in the top panel represent the traveling waves at each end of the path at t = 0. 

times [m/N = 3/2) to assume the shape of the other hump at f = 0. The transition from 
this behavior to the braided behavior occurs at the point on path e that a third hump 
becomes recognizable in the wave profile. The solutions on the paths /, /i, i, j, k and 
I are similar to those on path e, but the braiding patterns are more complicated near the 
right end-points of these paths. 

All the traveling waves we have described until now move left. To see what happens 
to a right-moving wave, we computed the first bifurcation from the simplest such case and 
obtained the path 

(1,1,1,2) ^ (2,0,1,2). (43) 

Thus, the one-hump right-traveling wave is connected to the two-hump stationary solution. 
Solutions near the left end of this path consist of a large-amplitude, right-moving soliton 
traveling over a small- amplitude, left-moving soliton. As we progress along the path, the 
amplitude of the left-moving soliton increases until the solitons cease to fully merge at 
t = T/4 and t = 3T/4. Instead, a dimple forms in the wave profile at these times and the 
solitons begin to bounce off each other, trading amplitude so the right-moving wave is larger 
than the left-moving wave. This type of behavior has also been observed by Leveque |22] for 
the KdV equation for solitons of nearly equal amplitude. Both types of behavior (merging 
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and bouncing off one another) are illustrated in Figure HI As we proceed further along this 
path, the solitons settle into a synchronized dancing motion without changing their shape 
or deviating far from their initial positions. Eventually the "dancing amplitude" becomes 
small and the non-trivial solution turns into a stationary two- hump solution. 

In order to guess a general formula for the relationship between two traveling waves that 
are connected by a path of non-trivial solutions, we generated two additional paths, namely 

(2,0,2,2) ^ (3,-1,1,2), 
(3,0,3,3) ^ (4,-1,1,3). 

After studying all the paths listed in ([38|) ~ (|H|) . we propose the following conjecture, which 
we prove as part of Theorem [3] in Section HI 

Conjecture 2 The four-parameter sheet of non-trivial solutions with bifurcation parame- 
ters {N ,v,n,m) coincides with the sheet with parameters (N' ,u' ,n' ,m') if and only if 

ifn<N: N' = N - n, v' 
ifn>N: N' = n + 1, u' 

By symmetry, we may interchange the primed and unprimed indices in either formula; 
thus, N' > N 44> n < N 4^ n' > N'. In most of our numerical calculations, A^' turned out 
to be larger than N. In the exact formulas of Section HI we find it more convenient to adopt 
the convention that N' < N since, in that case, all the solutions on the path connecting 
these traveling waves turn out to be A-particle solutions as described in Section 12.11 

Equations (|45p and (j46p are consistent with the bifurcation rules of Appendix [A] in that 

n<N, menu + NZ u' e Z, m' e {n' + + N"L, (47) 

n>N, me{n+l)u + N'L => u' e Z, m' e nv' + N'Z. (48) 

However, if the mean is held constant, they do not necessarily respect the requirements 
on ao listed in Figure [HI in Appendix Rl For example, if ao < 3, then (2,1,1,1) is a 
valid bifurcation, but the re-connection (1, 1, 1, 1) predicted by (|35|) is legal only if ao = 3. 
Interestingly, when we use our numerical method to follow the path of non-trivial solutions 
that bifurcates from (2, 1, 1, 1) with the mean ao = 1-2 held constant, it does not connect up 
with another traveling wave. Instead, as illustrated in Figure O as we vary the bifurcation 
parameter, the two humps (of the solutions labeled A,B,C) grow in amplitude and merge 
together until they become a single soliton traveling very rapidly on top of a small amplitude 
stationary hump. As the bifurcation parameter p = ai(0) approaches a critical value, the 
period T approaches zero and the solution blows up in L^(0, 27r) with the Fourier coefficients 
of any time-slice decaying more and more slowly. 

As another example, the bifurcation (3, 1, 1, 1) is valid when ao < 5 but the reconnection 
(2, 1, 2, 1) is only valid if ao = 5. If we hold ao < 5 constant, the solution blows up as we vary 
p = 02(0) from to a critical value. However, if we simultaneously vary the mean so that 
it approaches 5, we do indeed reach a traveling wave with bifurcation indices (2,1,2,1). 
To check this numerically, we started at (3,1,1,1) with ao = 4.8 (which has a = 



[N — n)v + m 

N ' 
(n + l)u — m 

N ' 



n' = N - I, m' = m, (45) 

n' = n + 1 — A^, m' = m. (46) 



16 




Figure 4: Periodic solutions with mean ao = 0.544375 between the one-hump right-travehng 
wave (dotted curve, top panel) and the two-hump stationary solution (dotted curve, bottom 
panel). Top: a large, right-traveling soliton temporarily merges with a small, left traveling 
soliton at t = ^ and t = |T. Bottom: two solitons traveling in opposite directions bounce 
off each other at j and |T and change direction. 




Figure 5: Left: path of non-trivial solutions with mean ag = 1-2 that bifurcates with 
indices (2,1,1,1) from the two-hump traveling wave. These solutions do not re-connect 
with another traveling wave, but instead blow up as T — t- 0. The solution Q is shown at 
right, where a large, right-moving soliton travels rapidly over a small, stationary hump. 
The dotted curves are initial conditions for the points labeled A, B, C at left. 
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I (2, 1,2,1), P = ^1 consists of a traveling wave inside a football-shaped envelope. The 
exact solution appears to be of the form it(x, t) = A + B (sin |) sin — "^t) . 



|/3| = l/\/3l) and computed 40 solutions varying p from to 0.1 and setting uq = 4.8 + 2p. 
The bifurcation at the other end turned out to be (2, 1, 2, 1) with = b, f3 = \p = 0.025, 
a = (1 — 3/3^)/(l — /3^), T = tt/{5 — 2a), as predicted by Conjecture [2l The solutions on 
this path have the interesting property that the envelope of the solution pinches off into 
a football shape at one point in the transition from the three-hump traveling wave to the 
two-hump traveling wave. Using a bracketing technique, we were able to find a solution 
such that the value of u{0, t) remained constant in time to 8 digits of accuracy. The result 
is shown in Figure [6l 

In summary, it appears that the family of bifurcations with indices (A^, z^, n, m) is always 
connected to the family with indices {N', u', n', m') given by (j45p and ()46p by a sheet of 
non-trivial solutions, but we often have to vary both the mean and a Fourier coefficient of 
the initial condition to achieve a re-connection. Thus, the manifold of non-trivial solutions 
is genuinely two-dimensional (or four dimensional if phase shifts are included) . Some of its 
important properties cannot be seen if we hold the mean uq constant. 

4 Exact Solutions 

In this section we use data fitting techniques to determine the analytic form of the numerical 
solutions of Section [3l We then state a theorem that confirms our numerical predictions 
and explains why some paths of solutions reconnect with traveling waves when the mean is 
held fixed while others lead to blow-up. The theorem is proved in Appendix [Bj 

4.1 Fourier CoelRcients and Lattice Sums 

One striking feature of the time-periodic solutions we have found numerically is that the 
trajectories of the Fourier modes Ck{t) are often circular or nearly circular. Other Fourier 
modes have more complicated trajectories resembling cartioids, fiowers and many other 
familiar "spirograph" patterns (see Figure [T]). This led us to experiment with data fitting 
to try to guess the analytic form of these solutions. The first thing we noticed was that the 
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trajectories of the spatial Fourier coefficients are band-limited in time, with the width of 
the band growing linearly with the wave number: 



u{x,t)= ^ c,.(^)e^'=^ Ck{t)= Ckje-'^^\ c^j = if |j| > r|A:|. (49) 

k=—oo j=—oo 

Here r is a fixed positive integer (depending on which path of non-trivial solutions u belongs 
to) and the c^j are real numbers when a suitable choice of spatial and temporal phase is 
made. Since u is real, these coefficients satisfy C-k-j = Ckj- 

Each path of non-trivial time-periodic solutions has a lattice pattern of non-zero Fourier 
coefficients Ckj associated with it. In Figure El we show the lattice of integers {k,j) such 
that Ckj / for solutions on the paths 

(1,0,1,1) ^ (2,-1,1,1), (2,-1,2,3) ^ (3,-3,1,3), 
(1,1,1,2) ^ (2,0,1,2), (2,-1,4,11) ^ (5,-8,3,11). 

All solutions on a given path have the same lattice pattern (of solid dots), but different 
paths have different patterns. One may show that if u{x, t) is of the form (I49p and 

k X / 27r \ 

2 Z_/'^'P'^'=-'J-P = ( ^1^1 + ) "^^^i' {k>0,jeZ), (51) 
i,p ^ ^ 

then u{x, t) satisfies the Benjamin-Ono equation, uu^ = Hu^x — ut- The traveling waves at 
each end of the path have fewer non-zero entries, namely 



Ckj 



( Airy 2:21k. I. ^ j 

^^"^ NT K J U, l-3/3^\ 

2iV/3lWA^ keNZ\{0}, j = !§ \, (^« = y3^]- (52) 



otherwise. 



Here a tilde is used to indicate a solution about which we linearize. Substitution of Ckj = 
Ckj + Edkj into (j5ip and matching terms of order e leads to an eigenvalue problem with 
solution 

'zM,n{k), kekN,n+NZ, i = ^, 

dkj = I ZN,n{-k), k G -kN,n + NZ, j = (53) 
^0 otherwise, 

with zis[,n{k) as in (jlSp . The non-zero coefficients dkj in this linearization are represented by 
open squares in Figure [HI Recall from (fTSj) that if n > and k < n — N then ZN,n{k) = 0, 
but iin< N, the non-zero entries oi zj^ ^ik) continue in both directions (with k approaching 
+00 or — oo). This is why the rows of open squares terminate in the graphs in the top row 
of Figure [8] rather than continuing past the origin as in the graphs in the bottom row. 



4.2 Elementary Symmetric Functions 

It is interesting that the lattice patterns that arise for the exact solutions (beyond the 
linearization) contain only positive integer combinations of the lattice points of the lin- 
earization and of the traveling wave (treating the left and right half-planes separately). 
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Figure 7: Left: Trajectories C2{t) for five solutions on path g in ()4ip . The evolution of C2{t) 
on paths / and h in ()4ip are similar, but with three- and five-fold symmetry rather than 
four. Right: Trajectories Cs{t) for three solutions on path i in (j4ip . 




Figure 8: Each pair (aligned vertically) corresponds to a path of non-trivial solutions 
connecting two traveling waves. Solid dots represent the non-zero entries c^j in ()49p of the 
exact solutions along this path; open circles represent a traveling wave; and open squares 
represent the non-zero entries dkj in the linearization about the traveling wave. 
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Somehow the double convolution in (|5ip leads to exact cancellation at all other lattice sites! 
This suggests that the Ckj have a highly regular structure that generalizes the simple power 
law decay rate of the Fourier coefficients nstat(^) -^j of the A^-hump stationary solution. 

The first step to understand this is to realize that there is a close connection between 
the trajectories of the Fourier coefficients and the trajectories of the elementary symmet- 
ric functions of the particles f3i, f^N in © above. Specifically, because the Fourier 
coefficients of 4'{x] j3) in (j4j) are of the form 2(3^ for /c > 1, we have 



Pl{t) + --- + P%{t) = -cu{t), 



k>l, Ckit) 



1 

2^ 



2n 



u{x, t)e 



-ikx 



dx 



Next we define the elementary symmetric functions cjj via 



(j = i,---,iv) 



(54) 



(55) 



so that 

It is well known 
/ 1 

E = 



\±aN ■ ■ 



P{z) 



N 

1=1 



N 



'aiZ 



N-j 



j=0 

that the companion matrix E of P has the Jordan canonical form 



(56) 



1 
-£72 ai I 



V 



Jm / 





V 



1 









1 



where I : {1, . . . , m} — )• {1, . . . , N} is an enumeration of the distinct roots of P{z) = and 
the size of the Jordan block Jr is equal to the multiplicity of As a result, the trace of 
powers of S will give the power sums of the /3i, and hence the Fourier coefficients: 



Cfc = 2tr(S'=), 



(A;>1). 



(57) 



Thus, if the elementary symmetric functions are finite sums of circular orbits, then the 
Fourier coefficients will be as well, and we expect higher Fourier modes to involve more 
terms, in accordance with our findings above. 

Before presenting our main result, note that once the mapping (I45p from (A^, i^, n, m) 
to {N' ,1/' ,n' ,m') is known, we can choose A^, u, N' and z^' independently, subject to the 
conditions 

N' 



N' < N, 



V > — u. 

N 



(58) 



The first condition is merely a labeling convention while the second is an actual restriction 
on which traveling waves are connected together by a path of non-trivial solutions. The 
formulas of Conjecture [2] then imply that 



m = m' = Nv' - N'u > 0, 



n 



N -N', 



n 



N-l. 



(59) 
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After extensive experimentation with data fitting on the numerical simulations described in 
Section [3l we arrived at the form (|6ip below for the polynomial P. We then substituted the 
ansatz ([5U|) into ([T]) to obtain algebraic relationships between ^4, B, C, oq, u, N, N', v and 
v' ^ namely (|9ip - ()93p in Appendix|Bl We solved these using Mathematica to obtain formulas 
for yl, B and io in terms of C, ao; -^'i ^ ^-iid v' . We had to break the analysis into three 
cases depending on whether v is less than, equal to, or greater than v' . By comparing our 
exact solutions with previously known representations of multi-periodic solutions [26], we 
found that all three cases could be unified by replacing C and by two new parameters, p 
and p', related to C and ao by (|62]) below. We give a direct proof of the following theorem 
in Appendix [Bl 

Theorem 3 Lei N, N' , v and v' he integers satisfying N > N' > and Nu' - N'v > 0. 
There is a four-parameter family of time-periodic solutions connecting the traveling wave 
bifurcations {N' , N — l,m) and {N,v,N — N' ,m), where m = Nv' — N'v. These solutions 
are of the form 



N 



n(x,t) = ao + ^0(x;/3/(t)), ^{hl3) 



1=1 



'2^l^l, /t<0, 
0, A: = 0, (60) 

2p\ k>0, 



where . . . , /3Ar(t) are the roots of the polynomial 

P{z) =z^ + Ae-^^'^h^-^' + Be-'^^-^'^^h^' + Ce"*'"^* (61) 

with 



iN-N'-^p + p' / {N + p')p' 



N + p + p' \i N'{N - N') + {N + p')p'' 



i{u—u')Luto —i{N—N')xQ 



N'{N - N') + {N + p')p'\/ N-N' + p' 



(62) 



iN-N' + p + p' 



N -N' -^py N + p + p' 



N'^u'-{N'fu 2N'{u'-v), 27r N' {N - N'){N + 2p') 

ao = 2p p, uj = — = . 

m m I m 

The four parameters are p > 0, p' > 0, xq GM. and to ^ The N- and N' -hump traveling 
waves occur when p' = and p = 0, respectively. When both are zero, we obtain the constant 
solution u(x,t) = I^LJijdJUji _ 

Remark 4 The parameters xq and to spatial and temporal phase shifts. A straightfor- 
ward calculation shows that if u has parameters p, p' , xq and to in Theorem [3] while u has 
parameters p, p' , and 0, then u{x,t) = u{x — XQ,t — to). 

There are two features of this theorem that are new. First, it had not previously 
been observed that the dynamics of the Fourier modes of multiperiodic solutions was so 
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simple. And second, in our representation, it is clear that these solutions reduce to traveling 
waves in the limit as p or p' approaches zero. By contrast, other representations become 
indeterminate in the equivalent limit, and are missing a key degree of freedom (the mean) 
to allow bifurcation between levels of the hierarchy of multi-periodic solutions. 



4.3 Three Types of Reconnection 

We now wish to explain why following a path of non-trivial solutions with the mean uq held 
fixed sometimes leads to re-connection with a different traveling wave and sometimes leads 
to blow-up of the initial condition. By Theorem [3l ao depends on the parameters p and p' 
via 

, 2N'{v' - v) , , N^v' - {N'fu 

ao = dn — 2p p , On := . (bo) 

m m 

If we hold fixed, then p and p' must satisfy 

^p^ 2N'{u'-u) ^,^ (64) 
m 

This is a line in the p-p' -plane whose intersection with the first quadrant gives the set of 
legal parameters for a time-periodic solution to exist. We assume the mean is chosen so 
that this intersection is non-empty. If the p- or /^'-intercept of this line is positive, the 
corresponding traveling wave bifurcation exists. There are three cases to consider. 

Case 1: [v < v') Both intercepts will be positive as long as oq < Og- Thus, a reconnec- 
tion occurs regardless of which side of the path we start on. 

Case 2: {u = v') The line (j64p is vertical in this case, so p = (og — ao)/2 remains 
constant as we vary p' from to oo. As — t- oo, we see from (|62|) that T — t- 0, vl — )• 1, and 
B and C both approach ^J pj (N — N' + p). In this limit, A^' of the roots /?/ lie on the unit 
circle at t = 0, indicating that the norm of the initial condition blows up as p' — t- oo. 

Case 3: {u > u') The line has positive slope in this case, li < a^, a, bifurcation 
from the A'''-hump traveling wave exists. If ckq > Oq, a bifurcation from the A^-hump 
traveling wave exists. And if oq = ckq) a bifurcation directly from the constant solution 
u = Oq to a non-trivial time periodic solution exists. In any of these cases, another traveling 
wave is not reached as we increase p and p' to oo. Instead, T — t- and A, B and C all 
approach 1. As a result, all the roots f3i approach the unit circle, indicating that the norm 
of the initial condition blows up as p, p' — )• oo. 

Example 5 Consider the three-particle solutions on the path e : (2, —1, 2, 3) -f-T- (3, —3, 1,3) 
in Figures [T] and [3l Since —3 = < i^' = — 1, we do not need to vary the mean in order to 
reconnect with a traveling wave on the other side of the path. Suppose < = lis held 
fixed. Then the parameters p and p' in Theorem [3] satisfy 

If 8 A , 3(1 - ao) 

p=-[l-ao--p'y ^<p'< ^ g °^ (65) 

The solutions u{x,t) on this path are of the form ([60]) with particles /3/(i) evolving as the 
roots of the polynomial 

P{z) = z^ + Ae'^^h + Be^'^^z'^ + Ce^'^\ (66) 
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particle trajectories particle trajectories particle trajectories particle trajectories 




Figure 9: Trajectories (3i{t) for four solutions on the path (2,-1,2,3) -f-)- (3,-3,1,3) with 
mean ao = 0.544375. The markers give the position of the /3/ at t = 0. The value of p' in 
dSSD is, from left to right: 0.1707, 0.1642, 0.1634 and 0.1369. In Figure El p' = 0.0862. 



where 

/ (9-3ao-2pO(3 + pOp^ ^ / (3 - 3ao - 8p')(3 + p')p' 



21 - 3ao - 2p'){2 + + p')' V - 3ao - 8p'){2 + + P') ' 

. (67) 

^_ / (9-3ao-2p0(3-3ao-8p0 ^ _ 2^ _ 2(34^V) 
~ y (21 - 3ao - 2p'){9 - 3ao - 8/)') ' ^ - ^ - 3 

The transition from the two- to three-hump traveling wave occurs as we decrease the bifur- 
cation parameter p' from 3(1 — ao)/8 to 0. This causes C to increase from to ^J and 



A to decrease from y i9_^3^q to 0. B is zero at both ends of the path. 

The trajectories /3i(t), /?2(i) and /33(t) for ao = 0.544375 and four choices of p' are shown 
in FigurelH For this value of the mean, p' varies from 0.17086 to 0. Note that the bifurcation 
from the two-hump traveling wave causes a new particle to nucleate at the origin. As p' 
decreases, the new particle's trajectory grows in amplitude until it joins up with the orbits 
of the outer particles. There is a critical value of p' at which the particles collide and the 
solution of the ODE ([3|) ceases to exist for all time; nevertheless, the representation of u in 
terms of P in (j83|) in Appendix [B] remains well-behaved and does satisfy dH) for all time. 
Thus, a change in topology of the orbits does not manifest itself as a singularity in the 
solution of the PDE. As p' decreases further, the three orbits become nearly circular and 
eventually coalesce into a single circular orbit (with i/ = — 3) at the three- hump traveling 
wave. The "braided" effect of the solution shown in Figure [3] is recognizable for p' < 0.15 
or so for this value of the mean. 



5 Interior Bifurcations 

We conclude this work by mentioning that our numerical method for following paths of 
non-trivial solutions from one traveling wave to another occasionally wanders off course. 
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following an interior bifurcation rather than reaching the traveling wave on the other side 
of the original path. These interior bifurcations lead to new paths of non-trivial solutions 
that are more complicated than those on the original path. For example, on the path 

(1,1,1,2) ^ (2,0,1,2), (68) 

Theorem [3] tells us that the exact solution is a two-particle solution with elementary sym- 
metric functions of the form 

aiit) = -{Ae~'^' + Be'^'), a2{t) = C. (69) 

We freeze ao < Oq = 2, set /? = ^(2 — ao — p')^ and determine that 



A = ,-^{x.-^t,) {^-ao + p'){2 + p')p' „ ^ -i(xo+-to) {2-ao-p')i2 + p')p' 



y (6-ao + /o )(4-ao-p') T 

In Figure [TUl we show the bifurcation diagram for the transition from the one-hump right- 
traveling wave (labeled P) to the two- hump stationary solution (labeled Q). This diagram 
was computed numerically before we had any idea that exact solutions for this problem 
exist; therefore, we used the real part of the first Fourier mode at t = for the bifurcation 
parameter rather than p' . We can obtain the same curves analytically as follows. The upper 
curve from P to Q (containing A1-A5) can be plotted parametrically by setting xq = 7r/2 
and to = vr/2a; in (f70]l . varying p' from 2 — to 0, holding ao = 0.544375 fixed, and plotting 
—2{A + B) versus T = jq^- The lower curve from P to Q is obtained in the same fashion 
if we instead set xq = to = 0. 

As illustrated in Figure [TOl solutions such as A1-A5 on the upper path have cTi(t) 
executing elliptical, clockwise orbits that start out circular at the one-hump traveling wave 
but become more eccentric and collapse to a point as we progress toward the two-hump 
stationary solution Q. Meanwhile, cT2(t) remains constant in time, nucleating from the 

origin at the one- hump traveling wave and terminating with a2 = —•^J^E^ at the two- 
hump stationary solution. On the lower path, the major axis of the orbit of cJi is horizontal 
rather than vertical and (T2 moves right rather than left as we move from P to Q. 

When computing these paths from P to Q, we encountered two interior bifurcations. 
In the bifurcation labeled B6 in Figure [TOl an additional elementary symmetric function 
nucleates at the origin and the trajectories of o"i and o"2 become more complicated. Through 
data fitting, we find that 

ai(t) = -(Ae-^"^* + Be''^^ + Cie^''^^), (71) 
a^it) = C + Cae^^"* + Cge^^'^*, (72) 
a-sit) = -C^e^"^\ (73) 

where the new coefficients Cj are all real parameters. We have not attempted to derive 
algebraic relationships among these parameters to obtain exact solutions. These trajectories 
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Figure 10: Left: Bifurcation diagram showing several interior bifurcations on the path 
(1,1,1,2) — )• (2,0,1,2). Right: Trajectories of the elementary symmetric functions cri(t), 
which have elhptical, clockwise orbits, and o"2(t), which remain stationary in time, for the 
solutions labeled A1-A5 in the bifurcation diagram. 




-0.2 0.2 -0.6 -0.4 -0.2 0.2 

Re(z) Re(z) 



Figure 11: Left: Trajectories of cri{t) for solutions labeled B1-B13 in Figure [TOj Center: 
Trajectories of a2{t) and a^lt). Since B6 is on the original path from P to Q, o"2(t) is 
constant and o'3(t) = for this solution. Right: The interior bifurcation causes additional 
lattice coefficients Ckj to become non-zero; grey circles represent the new terms. 
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are shown in Figure [TT] for the solutions labeled B1-B13 in the bifurcation diagram. The 
additional term in (j7ip causes the elliptical orbit of C7i(i) to deform by bulging out in the 
vertical and horizontal directions while pulling in along the diagonal directions (or vice 
versa, depending on which direction we follow the bifurcation). Meanwhile, cr2{t) ceases 
to be constant and (T^{t) ceases to be zero. To avoid clutter, we plotted the trajectories 
(72 (i) and (T3(t) for B1-B6 separately from B6-B13, illustrating the effect of following the 
bifurcation in one direction or the other. The additional terms in (I7ip - (j73p cause the lattice 
pattern of non-zero entries Ckj = ^ /q^ Cfc(t)e*-'^* to become more complicated, where we 
recall that in this case, 



ci^(t) = —[ u(x,t)e-^'=^dx = 2tr 
271" Jo 



10 
1 

y<J3{t) -CT2(t) ^Jl(t), 



The solid dots in Figure [TT] represent the non-zero entries of solutions on the original path 
from P to Q while grey circles show the additional terms that are non-zero after the bi- 
furcation at B6. Although this bifurcation causes some of the unoccupied lattice sites to 
be filled in, the new lattice pattern is rather similar to the original pattern and maintains 
its checkerboard structure. Also, this bifurcation leads to symmetric perturbations of the 
Fourier mode trajectories, and is also present (in a phase shifted form) along the lower path 
from P to Q. 

In the bifurcation labeled CO in Figure [TOt the fill-in pattern of the lattice representation 
is much more complicated, and in fact the checkerboard structure of the non-zero coefficients 
Ckj is destroyed; see Figure [T2j But actually, the elementary symmetric functions behave 
similarly to the previous case: By fitting our numerical data, we find that 

ai(t) = -{Ae-''^' + Be'^' + Cie^'^'), (74) 
a2 {t) = C + C72e3^-* + Cse^'"* , (75) 
(73 (t) = -Qe^^"*, (76) 

so each of the new terms executes one additional loop per cycle of the periodic solution in 
comparison to the corresponding term in (ITT!) ~ ([73]) . This extra loop causes a star-shaped 
perturbation of the ai ellipse instead of the rectangular and diamond shaped perturbations 
seen previously in Figure [TTJ As a result, this bifurcation is not present on the upper path 
from P to Q because the symmetry of the perturbation does not respect the 90 degree 
rotation of the orbit ai{t) associated with the ^-spatial and ^-temporal phase shifts that 
relate solutions on the upper and lower paths from P to Q. To follow the bifurcation at CO 
in the other direction, we can use the same numerical values for A, B, C, Ci, C2, C3, C4 in 
(|74 p - (j76p after changing the signs of the latter four parameters. This causes the trajectories 
of cji in Figure [12] to be rotated 180° with a corresponding ^ phase-shift in time so that 
the initial position o"i(0) remains on the left side of the figure. Meanwhile, the trajectory 
of £72 (i) experiences a ^ phase-shift in time with no change in the location of the orbit, and 
(73 (t) starts on the opposite side of its circular trajectory about the origin. 

In Figure [T3] we show the orbits of the 16th and 26th Fourier modes for the solution 
labeled C9 in the bifurcation diagram of Figure [10] As the index of the Fourier mode in- 
creases, these trajectories become increasingly complicated (involving more non-zero terms 
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-0.75 -0.5 -0.25 0.25 0.5 

Re(z) 

Figure 12: Left: This interior bifurcation causes more lattice coefficients to become non- 
zero than the interior bifurcation of Figure [HI Right: Trajectories of cri(t), (72 (t), and cr^^t) 
for the solutions labeled C0-C9 in Figure [TU) The long axis of the ellipse CO is horizontal 
because we start from the bottom branch connecting P to Q in Figure [TOl 




Figure 13: The trajectories of the Fourier modes become very complicated after the interior 
bifurcation occurs. Here we show the 16th and 26th Fourier modes Cfc(t) = ak{t) + ibk{t) 
over one period. It was clearly essential to use a high order (in fact spectrally accurate) 
numerical method to resolve these dynamics when computing time-periodic solutions. 
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Ckj in the lattice representation), but also decay exponentially so that the amplitude of the 
orbit is eventually smaller than can be resolved using floating point arithmetic. We empha- 
size that these trajectories were resolved to full machine precision by our general purpose 
numerical method for finding periodic solutions of non-linear PDE (without any knowledge 
of the solitonic structure of the solutions). Everything we learned about the form of the 
exact solutions came about from studying these numerical solutions, which was possible 
only because our numerical results are correct to 10-15 digits of accuracy. 



A Bifurcation formulas and rules 



In this section we collect formulas relating the period, mean and decay parameter at a 
bifurcation. We also identify bifurcation rules governing the legal values of ao foi' ^ given 
set of bifurcation indices. 

In computing the nullspace J\f = kev DF{Uo,T) in Section [2131 we considered A^, z^, /3, 
T (and hence oq) to be given and searched for compatible indices n and m. The decay 
parameter |/3|, the mean oq, and the period T cannot be specified independently; any two 
of them determines the third. We now derive formulas for the period and mean in terms 
of (A^, I/, n, m) and /3. To simplify the formulas, we work with a = (1 — 3|/3p)/(l — |/3p) 
instead of /?. Note that as we increase |/3| from to 1, a decreases from 1 to — oo. For the 
period, we have 

27rm 

2vrm 

27rm 



n < N, 
n> N, 



(77) 



I N{n + 1- N){n + 1 + N{1- a)) 

so the period is independent of f3 when n < N, and otherwise decreases to zero as 
from to 1. For the mean, ao, we note that 

27ri/ 



cT 



Hence, using 



N ' 

— we obtain 



Na 



ao = Na + 



27rz^ 
iVT' 



varies 



(78) 



n(N — n) 
N + — -u 



N + 



m 

(n + 1 



(l-a)iV, 



N){n + 1) 



m 



1 



n+l- N 



m 



N{1 



a 



n< N, 



n> N. 



(79) 



Thus, as 1/3 1 varies from to 1, the mean ao decreases to — oo if n < N , and otherwise 
either decreases to — oo, increases to -|-oo, or is independent of /3, depending on the sign of 
[m-{n+l- N)v]. 

In practice, we often wish to start with A^, i/, n, m and ao and determine T and |/3| 
from these. However, not all values of ao are compatible with a given set of indices. The 
bifurcation rules are summarized in Figure [141 Solving (|79p for a yields 



a 



1 



(A^ — ao)?7i + n{N — n)v 
Nm ' 
(A^ - ao)m + {n + l- N){n + l)u 

[m-{n+l- N)u]N ' 



n< N, 
n> N. 



(80) 
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1. iV > 1, G Z, n > 1, m > 1 

2. if n < iV then 

• m £ nv + NIL 

n(N — n) 

• aQ<N + — -V 

m 

3. if n > then 

• m£ {n+l)iy + NZ 

• \im> {n + l- N)v then < N + ^ 

• if m < (n + 1 - N)v then ao > + i 

• if m = (n + 1 — A^)zv then = n + \ + N 



Figure 14: Bifurcation rules governing which values of uq are compatible with the bifurcation 
indices {N^v,n,m). 



The corresponding period is given by 

( 27rm 



T 



Nn{N - n) 
27r 



n+l-N 



I N{n + 1 + N -ao)'' 



n< N, 

(81) 

n> N. 



In the indeterminate cases {n > N, m = (n + 1 — A^)i^, = ^t- + 1 + N}, any a < 1 is 
allowed and formula (j77p should be used to determine T. 

If we express n, n', m and m' in terms of N, N' , ly' , then (|77p and (|79p give 

27r(jV^^ - iVV) , , Ar2^/ _ (iv/)2^ 

N'{N - N')N ' u V ; , ^ Nu' - N'u 

= N'{N-N')[N + {l-a')N'y "° = "° " iV.- - iVV ^ " " ^' ^''^ 

where q = and a' = ^l^j^/p ■ We note that the two traveling waves reduce to the 

same constant function when /3 — )• and /?' — )• 0, which is further evidence that a single 
sheet of non-trivial solutions connects these two families of traveling waves. 

B Proof of Theorem [3] 

As explained in Remark HI xq and to are spatial and temporal phase shifts, so we may set 
them to zero without loss of generality. We can express the solution directly in terms of the 
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elementary symmetric functions via 
N 



N 



ikx 



zd,P{z) 



(83) 



1=1 

N 



1=1 



.k=l 



+ ^4Re 



1=1 



1} = ao + 4 Re 



P{z) 



N 



Next we derive algebraic expressions relating A, B, C, ao, co, N, N', v and v' by substituting 
(|83p into the Benjamin-Ono equation ([1]). To this end, we include the time dependence of 
P in the notation and write (1831) in the form 



+2 



-idrh 



N 



where 



g{x,t) = P{e-'-,e 



ix —iu)t\ 



h{x,t) = g{x,t), 



P{z, \) = z^ + AX"' z"-''' + BX"'"' z''' + CX\ 



,u' N-N' 



v~v' N' 



(84) 



(85) 
(86) 



Note that P is a polynomial in z and a Laurent polynomial in A (as and v' may be 
negative). We may assume a; > 0; if not, we can change the sign of w without changing 
the solution by replacing (A, B, v, v', N') by {B, A, —u, u' — v,N — N'). Assuming the roots 
/3i{t) of z I— 7- P{z, e"*"^*) remain inside the unit disk A for all t, we have 



idx9 



N 



N oo 



ikx 



Hu 



1=1 k=l 



dx9 



+ Ni\+2 



Ni 



(87) 



Using (|84p and "9*^^^^ = ^^(^)' technique we learned by studying the bilinear for- 
malism approach of [32l |26]), the equation ^ {ut — Huxx + uu^) = becomes 



dx 



9 h 



dxQ dxh\ 1 ( , ^. ( dxQ dxh 

— + ^ + 7 «o - 4Af + 2« — - ^ 
9 h J 4\ \ g h 



0. 



The expression in brackets must be a constant, which we denote by 7. We now write 

p,k = izd,y{xdxfPiz,x) 



z=e 



SO that e.g. dt9 = —iojPoi and dxh = iPiQ. Equation (I88p then becomes 

T^ooi^oo + i^oo [^^20 + ujPqi + (ao - 4iV)Pio] 

+ Poo [P20 + ooPoi + (ao - 4iV)Pio] + 2PioPio = 0, 

where we have absorbed ^(ao — 4A^)^ into 7. This equation may be written 



(90) 



+ e4z'' X-" +65 = 0, 
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where [aj = a + a = 2 Re{a}, 

ei = [-i + vuj + N"^ + (ao - 4iV) A^] C, 62 = [7 + J^w + iV^ + (ao - 4A^) A^] AB, 
and, after setting 7 = (3A^ — aQ)N — voj to achieve ei = 62 = 0, 

63 = [(Ar')2 _ 2NN' + A^'ao - v'^]B + [{N'f + 2iVAf' - N'uq + z^'wJAC = 0, (91) 

64 = [SiV^ - 4A^iV' + {N'f -{N- N')ao + - u')u;]BC 

- [N^ - {N'f -{N - N')ao + {u - iy')uj] A = 0, (92) 

65 = {Nao -uoj- iV^) + [{2N' - N)ao + {u - 2v')uj + SiV^ - SNN' + 4(iV')^]-B2 

+ [{N - 2N')ao + 4.{N'f - + {2v' - u)uj\ + [(3A^ - ao)iV + vu] = 0. (93) 

Using a computer algebra system, it is easy to check that (|9ip - ()93p hold when A, B, C, ao 
and u! are defined as in ( 1621) . When n' = 0, we have A = B = and C = . — - — so that 

y N + p 



2 



R(f\ ^r-rV- ^r-r.-^ct , N'{N-N')v ^ A-'iC 

Pi[t) = V-CA'^= V-Ce , c= — = = ao-N- 

iV m 1 — 

where each /3; is assigned a distinct A^th root of — C. By ([5]), this is an A^-hump traveling 
wave with speed index u and period T = —. Similarly, when p = 0, we have B = C = 

and A= \ \ — so that 



V=^e-*^* 1<N'\ ujv' {N-N'){N + 2p'y ,wl-3^' 
Piit) = < n j^^T>}' c= — = = ao-N- 



/ > AT' j ' A^' m " 1 - ^2 ' 

which is an A^'-hump traveling wave with speed index u' and period T = 

Finally, we show that the roots of P(-, A) are inside the unit disk for any A on the unit 
circle, S^. We will use Rouche's theorem [Tj. Let 

/,(^) = z^+^A'^'z^-^' + sy-'^'z'"' + CA^ 

/2(.) = z^+^A'^'.^-^', 

f3{z) = z^ +By-'''z^'. 

From ([62]), we see that {A, B, C} C [0, 1), A > BC, B> CA and C> AB. Thus, 

d2{z) := |/2(z)|2 - \Mz) - /2(z)P = lA-'^'z^' + A\^ - \BX-'z^' + 

= 1 + ^2 _ ^2 _ ^ 2(^ - 5C) cos > (1 - A)2 - (B - C)^ 

where \~^' z^' = e*^. Similarly, 

ds := |/3(z)P - \h{z) - Hz)? > (1 - Bf -{A- Cf. (95) 
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Note that 

B<A, C<B B-C<B-AB<1-A ^2(2) > for z G 5\ 

B<A, OB => \C-A\<l-B d3(z) > for z G 5\ 

A<B, C<A A-C <A- AB <l- B dsiz) > for z e , 

A<B, OA =^ \C-B\<1-A d2{z) >OfoT z e S^. 

Thus, in all cases, fi{z) = P{z, A) has the same number of zeros inside as f2{z) or /^{z), 
which each have A'' roots inside S^. Since fi{z) is a polynomial of degree N, all the roots 
are inside S^. 
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